Heuristic derivation of continuum kinetic equations from microscopic dynamics 
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We present an approximate and heuristic scheme for the derivation of continuum kinetic equations 
from microscopic dynamics for stochastic, interacting systems. The method consists of a mean- 
field type, decoupled approximation of the master equation followed by the 'naive' continuum limit. 
The Ising model and driven diffusive systems are used as illustrations. The equations derived are 
in agreement with other approaches, and consequences of the microscopic dependences of coarse- 
grained parameters compare favorably with exact or high-temperature expansions. The method is 
valuable when more systematic and rigorous approaches fail, and when microscopic inputs in the 
continuum theory are desirable. 

PACS numbers: 05.50.+q, 05.10.Gg, 64.60.Cn 



I. INTRODUCTION 

Ever since its introduction in a classic treatment of the 
Brownian motion 111, the Langevin equation has been 
playing an important role in modern statistical physics. 
It provides a mathematical framework and a physical ba- 
sis for studying stochastic processes in statistical, me- 
chanical systems. Applications are wide-ranging Mi, in- 
cluding chemical reactions, laser physics, diffusive pro- 
cesses, and modern theories of dynamical critical phe- 
nomena (3) . Recent topics such as surface growth jj] and 
pattern formation j^] also rely heavily on the Langevin 
equation. 

It is fair to say, however, that despite its popularity, 
the Langevin equation for a specific problem is seldom 
derived from the corresponding microscopies. It is often 
postulated on grounds of symmetry and physical reason- 
ing. Only rarely in simple circumstances is it derived, for 
example, in reasonable details from the more fundamen- 
tal master equation. In this article, we shall present an 
approximate scheme for deriving the Langevin equation, 
starting from the microscopic specification of the dynam- 
ics. Our goal is not to provide a formal derivation, but 
rather to propose a simple and heuristic means that can 
be applied generally to many stochastic systems, several 
of which shall be discussed below. 

Expositions of the historical, philosophical and tech- 
nical aspects of the Langevin equation are beyond the 
scope of this article, interested readers are thus referred 
to the relevant literature 0. This paper is organized 
as follows: An elementary recapitulation of the master 

tuation and Langevin equation is presented in Section 
Section III contains several examples as illustrations 
of the method, as well as assessment of the quality of the 
approximation involved. Conclusion is given in Section 



IV. A discussion of the noise correlation is given in the 



Appendix. 



II. MASTER EQUATION AND 
TIME-DEPENDENT GINZBURG-LANDAU 
EQUATION 

In statistical physics, one of the most important ap- 
plications of the Langevin equation is in the theories of 
dynamical critical phenomena ||. Therefore, our discus- 
sion shall be cast in that language, although it should be 
obvious that the method itself is not limited to systems 
exhibiting those phenomena. For concreteness, consider 
the kinetic Ising model that obeys Glauber (i.e., spin- 
flip) dynamics At the classical, microscopic level of 
description, the system consists of N spins <r, interacting 
via the Hamiltonian 



H = - Jy <7 l (jj, 

(i,j) 



(1) 



where J is the coupling constant, and denotes a 

sum over nearest-neighbor pairs. The time evolution of 
the system is governed by a master equation 

P{a;t + l)-P{B;t) = 

M<? - 3)p&; - -> v')p& *)] . ( 2 ) 

a' 

where P{3; t) is the joint probability of finding the sys- 
tem in the spin configuration a = {<7i, 02, • • • , <7jv} a t 
time t, and w's are the transition rates between two con- 
figurations that differ only by one spin flip. There is a 
great deal of freedom in the choice of w, as long as the 
following detailed balance condition is satisfied to ensure 
the same equilibrium distribution P eq (a) ~ e^^ H ^ 



w(a' -> a) _ P cq (op 
w(a -> a') ~ P eq (a>) 



-0{H(a)-H{9>)] 



(3) 



with p = 1/ksT. In practice, the choice is largely dic- 
tated by mathematical convenience. The most common 
choices of W are the Metropolis rate in Monte Carlo sim- 
ulations for its ease of implementation, and the heat-bath 



1 



rate (also known as the Kawasaki rate 0) in analytic cal- 
culations for its analyticity. In this paper, we confine our 
attention to the latter choice. It is given by 



w(a — > a) 



1 



I + e -(3[H(3)-H{S')] ■ 



(4) 



Since the master equation is not very convenient for 
analytic purposes such as a renormalization group anal- 
ysis, one often turns to a mesoscopic, continuum repre- 
sentation. For an Ising system with Glauber dynamics, 
the relevant continuum field is the local magnetization 
density 6(f, t), which obeys a kinetic equation 



dl = _ r 6H 

dt 56 + ^' 



1 



H = I d d r \ -(V6) 2 + V(6) 



(5) 

(6) 
(7) 



This is an example of the time-dependent Ginzburg- 
Landau (TDGL) kinetic equation. In (||), d is the di- 
mensionality of the system, and TL is a coarse-grained 
Hamiltonian. For (||) to describe a stochastic process, 
the noise term £(r, t) is needed, which accounts for the 
effect of thermal fluctuations and prevents the system 
from trapping in metastable states. For mathematical 
convenience, it is often taken to be Gaussian with zero 



(C(r,t)) 

mm?,*)) 



2D5(r 



?)5(t-t'). 



(8) 
(9) 



For equilibrium systems, the correlation D in (^) has to 
be chosen to ensure that the stationary solution of Eq. (|5|) 



is consistent with the Boltzmann weight V, 



eq 



-H 



Appendix). For a system as simple as the Ising model, 
the static continuum Hamiltonian H can actually be de- 
rived from the microscopic H via the partition function 
by means of the Hubbard-Stratonovich transformation 
||, which is a trick based on the Gaussian integral. For 
the dynamics, the TDGL equation can be derived by 
coarse graining the master equation H|l(| , which in prin- 
ciple yields expressions of the coarse-grained parameters 
r, u and g in (@)-(R) as functions of microscopic ones in 

However, for more complicated H, the Hubbard- 
Stratonovich transformation fails and the coarse graining 
cannot be done explicitly. Moreover, these methods rely 
on the existence of a Hamiltonian H and the associated 
equilibrium Boltzmann weight e~^ H , which is not valid 
in generic non-equilibrium situations defined only by the 
dynamics jll| . For these reasons, it is highly desirable to 
have a way with which a continuum description can be 
extracted directly from the dynamics flI2| . 



III. FACTORIZATION AND NAIVE 
CONTINUUM EXPANSION 

Our method is very simple. It consists of two steps: 
a mean-field type factorization of joint probabilities into 
singlet ones in the master equation, followed by a 'naive' 
continuum expansion. The result is a continuum kinetic 
equation with full knowledge of the microscopic depen- 
dence of the coarse-grained parameters. Since the input 
is the master equation, whether the system is an equilib- 
rium one or not |fl| ] is irrelevant. To illustrate, we now 
discuss several examples in increasing order of sophisti- 
cation. 



A. ID Ising model 

Focusing on a spin at position i in a one-dimensional 
(ID) Ising model, it is easy to find by integrating out all 
other spins in Eq. (||) that 

P+(x; t+1)- P+(x; t) = P. w + P__ + w + 

+P+-+W+-+ + P+—W+— - P-+-W-+- 



P-A 



P++-W++- - P+++W+++, (10) 



where P+{x\t) denotes the singlet probability of finding 
the spin up at site x at time t, and P+++(ir; t) denotes the 
joint probability of finding three spins up at site x — 1, 
x, x + 1 respectively, and so on. From (|j), the heat- 
bath transition rates are given by w = = W4, 

w |_ = = w ^ + = = Wo, Wh h = 

w 1 = W-4, where 



W n 



1 



n(5J ' 



(ii) 



(cf. Adopting a mean-field approximation, the joint prob- 



abilities are replaced by their factorizations, e.g., 
P ++ -(x;t) P + (x - l;t)P + (x;t)P-(x + l;i). Since 
£ CT <7P CT (a:;t) = P+(x;t) - P-(x;t) = (a), and 

J2a P<y( x \t) = P+( x >t) + P-( x ',t) — 1) m the spirit of 
coarse graining we proceed to make the identification 



P±(x;t) 



(12) 



where <f> is the local magnetization density. By using 
4> instead of spin number densities, we take advantage 
of symmetries anticipated in the final kinetic equation. 
Since a spin flip depends on a total of z + 1 spins in (|^) , 
where z — 2d for hyper-cubic lattices, the factorization 
effectively produces a power series expansion in 6 up to 
6 Z+1 . After replacing P's by 0's, we make the transition 
to the continuum by 'naively' expanding about x, such 



(x±l;t) -> A(x;t) ± 



d6(x;t) ld 2 S(x;t) 



dx 



dx 2 



(13) 
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For most applications, we are only interested in the long- 
distance behavior, hence it suffices to stop at the lowest 
derivatives as shown. This procedure results in a deter- 
ministic kinetic equation for <fi in precisely the form of 
@"(0)j barring the noise term £: 



where the coefficients are given by 



1 



T=-{W-4-Wi), 



expanded in (3 J. Hence, our method has the advan- 
tage that it embodies a refined mean-field theory, 
as already applied to studies of stochastic resonance 
in Ising systems H]. 

6. Finally, due to the factorizations only the deter- 
ministic terms in (|^) can be derived. The noise 
term has to be deduced separately (see Appendix) . 
The result for the noise correlation D in (^) is 
D = k B T^T. 

Having gone through the details of our method, we 
now turn to a few less trivial examples. 



r = — (3W 4 - W-4 + 2W ) , 



9 = r= {Wi + W-i - 2W ) . 



(14) 



Several remarks are in order: 



1. Symmetries in the resulting continuum equation 
(with respect to <j), x and t) are as expected, be- 
cause the approximations respect those symmetries 
and leave them intact. 

2. There are explicit temperature dependences in the 
coefficients which cannot be deduced by symmetry 
or physical reasoning. Such dependences are spe- 
cific to the choice of jump rates which manifests 
through the approximations used. 

3. Noting that W- n = 1 — W„ for any n, we find T = 
\-W A > and g = at any T @, and r = 2W 4 /T 
has one zero, at T = 0. This is consistent with the 
absence of phase transition in the ID Ising model 
at any finite temperature, an improvement over the 
usual mean-field result T^ F = 2J/k B . There is no 
stability problem arising from g — because the 
quadratic coefficient is positive for T > 0. 

4. In the presence of an external magnetic field h, the 
degeneracies in jump rates are lifted (e.g., W+±_ = 
(1 + e ±2/3 ' l )~ 1 ). To 0(h), the kinetic equation ac- 
quires a new term Tfj,h on the right hand side, 
where 



23 

H = -fiWg + WiW-i). 



(15) 



Linear response then determines that the suscep- 
tibility is x = A*/ r = /3(1 — 7 2 /2)/(l — 7), where 
7 = tanh 2/3 J. In the Appendix we show that /1 is 
needed to fix the noise correlation. 

Besides capturing the correct symmetries, our re- 
sults compare quite well with exact results. From 
(||), the relaxation time can be read off easily as 
t = 1/Tr = 1/(1 — 7). This turns out to be exact 
H . For the susceptibility, deviation from the exact 
result !3e 2 P J [0 shows up only at 0((/3J) 5 ) when 



B. 2D Ising model 

The same procedure can be applied to the 2D Ising 
model with Glauber dynamics. Again we obtain (JH|), with 
the parameters given by 



T = -(-2W 4 + 2VF_ 4 - w 8 + Ws), 

8 

r = ^=(6W + 12W 4 - 4W- 4 + 5W$ - 3VF_ 8 ), 

8r 

g = -^-{~6W ~ 4W 4 + 4W_ 4 + 5W$ + W_ 8 ), 
H = ^(3VF 2 + 4WW-4 + WsWs)- 



(16) 
(17) 
(18) 
(19) 



It is worth noting that T, g and fi are positive definite for 
all T > 0, whereas r has one zero at T C GL 3.0898 J/k B w 
1.3616T C , again an improvement over the mean-field pre- 
diction T C MF = 4J/fc B , where T c = -2J/k B ]n(y/2- 1) « 
2.2692 J I k B is the exact critical temperature. As ex- 
pected, there is no <fi 5 and higher order term fl3[| . 

The results of r and % for the Gaussian case (g = 0) are 
quite satisfactory. They differ from high-temperature se- 
ries expansions [jl5| at order 0(((3J) b ) and 0(((3J) 4 ), re- 
spectively, whereas the usual mean-field results are worse, 
at 0{{l3Jf) and 0{{f3J) 2 ). 



C. 3D Ising model 

Despite being more tedious (128 terms on the right- 
hand side of the master equation), we also derive the 
kinetic equation for the 3D Ising model with Glauber 
dynamics. The results are: 



r = 



4M4 - 



451F4 



4W_ 8 + 5T^_ 4 - 5H/ 4 - 
W12), (20) 
18W_ 8 - 15H7 4 + 20Wo + 
30Ws + 7VFi 2 ), (21) 

+ 6VF_ 8 + 9VF_ 4 - 12W - 
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15W 4 + 6W 8 + 7W 12 ), (22) 
10W 2 + 15W_ 4 W 4 + 61^-8^8 + 
W- 12 W 12 ), (23) 



As for 2D Ising, T, g and fj, are positive definite for all 
T > 0, and r has one zero at Tp L « 5.0733 J/k B , 12% 
higher than the best estimate [fL6| , compared to T C MF = 
6J/fc B . 



D. ID driven lattice gas 

In many generic non-equilibrium systems [ ]Tl"| , the free 
energy does not exist and one has to start from the dy- 
namics, such as described by the master equation. A no- 
table example is the driven diffusive system |l7j , which is 
regarded as a paradigm of spatially extensive interacting 
systems that exhibit cooperative phenomena in steady- 
state non-equilibrium situations. In its standard form, 
it models an Ising-like lattice gas of particles whose mo- 
tion along a certain direction is biased by an external 
drive denoted by E. For E = 0, the model reduces to 
the ordinary kinetic Ising model with Kawasaki, or spin- 
exchange, dynamics (model B in ||). 

A question subject to recent debate concerns the form 
of nonlinearities associated with E (ljjl^] . That is an im- 
portant issue because the nonlinearities decide to which 
universality class of critical behavior the system belongs. 
It is interesting to see what the present method says 
about that. First, we consider a one-dimensional, sim- 
plified version in which the particles are not interact- 
ing except being hard-core, but their hoppings to nearest 
neighbors are biased by having different jump rates, p 
and q, to the right and left respectively. Hence, the mas- 
ter equation reads 

P + (x;t + l)-P+(x;t) = 

P P+-(x-l;t)+qP-+(x;t) 
- P P + -(x;t)-qP_ + (x-l;t), (24) 

where as usual an up(down) spin corresponds to the occu- 
pation of a particlc(holc), and joint probabilities such as 
P_i (x — 1; t) means the probability of finding a particle- 
hole pair at site x — 1 and x. After factorizations and 
applications of (12|) and (fi3f), we readily find 



d A = v d 2± 
at dx 2 



Sdf_ 
2 dx 



(25) 



where the diffusion coefficient is T> = (p + q)/2, as ex- 
pected, and the coefficient of driving is £ = (p — q). 
The nonlinear term is the same as in the 'standard' field 
theoretic model |20[| which was proposed on grounds of 
symmetries. As side remarks, note that we obtain the 
diffusion equation for p — q, and that £ is smooth in the 
'infinite' drive limit (p = 1, q = 0) which is used in most 
Monte Carlo simulations of driven diffusive systems. 



E. 2D driven lattice gas 

Generalization of the previous result to the 2D inter- 
acting driven lattice gas is immediate, despite the un- 
pleasant fact that there are altogether 512 terms in the 
master equation. In the presence of a drive E along the 
+y direction and attractive (J > in (|])) interaction 
between particles, the heat-bath rates for hoppings of 
particles along and against the drive take the form 



W, 



1 



n,±E 



I _|_ e n0J^ E j3J ' 



(26) 



where the dimensionless E (0 < E < oo) represents the 
"work done" on the particle by the field. Obviously, the 
rates for hoppings perpendicular to E are W n fi = W n . 

Going through the same procedure as above, we even- 
tually obtain a kinetic equation which is in complete 
agreement with the standard field theory of the driven 
diffusive system pi: 



5* 

dt ~ \ ax dx± 



dx 2 
£df_ 
2 W 



xy dx 2 dy 2 
d 2 \ I 



v d y 2 



y dy* 
d 2 



6 \ 9x dx 2 



<hj 



dy 2 



+ 



(27) 



The anisotropics are generated by the drive. Excluding 
the last term, this is the anisotropic generalization of the 
deterministic TDGL equation with conserved magnetiza- 
tion, i.e., model B S: 



dt 



(28) 



All coefficients are determined: 



a x = — - (69 - 85 W 4 - 68 W 8 - 17W 12 ) , (29) 

a xy = (20 -20W 4 - 16 W 8 -AW 12 + W- l2 .- E + 
256 

12, E + 4 VF_8,-£ + 4 W-s,E + 
5 W-4-E + 5 W-4.E - 5 W 4 -E - 5 W 4 ,e - 
AW^E-tWz,E-W l2 ,-E-Wi2,E), (30) 

ay = -^(8 W-v»,-b + 8 W- X2 ,e + 31 W- 8 ,-e + 
31 W- S ,e + 35 W- A ,- E + 35 W^e - 
10 Wq-e - 10 W ,e - 50 W^-e - 50 Wi, E - 
37^8,-£ - 37 W 8 , B - 9 W l2 ,- E - 9 W l2 . E ), (31) 

r x = ^(- 9 + 25 W i + 20 + 5 W 12 ), (32) 

r, = ^-(-2W_ 12 .- E -2W_ 12 , E -7W^ E - 
64 

7 W-$, E - 5 W-4,-s - 5 W- 4 , E + 10 Wq-e + 
10 W 0tE + 20 W4-E + 20 W 4 , E + 13 W S - E + 
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13Ws,e + 3W 12 ,-e + 3W 12 ,e), (33) 
9x = 4( 3 + w * ~ 4 w * ~ 3 Wia). ( 34 ) 

9y = 3^(6VK_i2,- J E + 6V^_i2, jE + 7^_ 8 ,_ £ ; + 
7 VF_ 8 , B - W_4,-B - VF_ 4 ,£ + 6 Wo.-E + 

6 Wo,e + 4 W 4 ,_e + 4 W 4 ,e ~ 13 Wg.-E - 
13W 8 ,e-9W 1 2,-e-9W 12 ,e) (35) 

f = "777 ( W—\2. — E + W-X2.E ~ 3 W- 8 ,-E + 

lb 

3 VF_ 8 ,B - 3 W-A-E + 3 W- AtE - 

2 Wo,-b + 2 VK , B - 3 W 4) -e + 3 W^ E - 

3Wg^E + 3W St E-W 12 ,-E + W 1 2 ! E). (36) 

They have these important properties: 

1. All but £ are even in E, consistent with the invari- 
ance of the dynamics under {E — > —E, y — ► — y}. 

2. The quadratic coefficient r x is independent of E. 
It has one zero at T C GL = 3.86143J/fcs, as shown 
in Fig. |l|. In contrast, r y depends on E. Fig. |2| 
displays the behavior of r y (T, E) versus r x (T) as T 
is lowered from above T C GL at fixed E, as well as 
r y (T = T C GL , E) versus E. It shows that for any 
E > 0, r x always vanishes before r y does when T 
is decreased. For small E, r y » r x + ci? 2 where 
c > 0. Consequently, at the critical temperature, 
the dominant derivatives come from the r y and a x 
terms, leading to the identification of an intrinsi- 
cally anisotropic critical theory with scaling of mo- 
menta k y ~ k x . This agrees with a previous per- 



turbativc argument 1 20 1 



3. The coefficient £ of the leading nonlinearity in- 
duced by the drive vanishes linearly in E at small 
E, and saturates to a constant at E = oo. This de- 
pendence, exhibited in Fig. |3|, is already anticipated 
above from the ID model and argued previously 
[ fL9[ , but at odds with the claims in . 



4. At E = 0: we find a x 



7^ cn xy - the continuum 



model derived is not rotationally invariant. This is 
not surprising because the initial lattice model is 
not. However, it turns out that a xy = a x at T GL . 
At present, we are not sure whether this acquire- 
ment of higher symmetry at the critical point is 
general for the underlying lattice model or specific 
to the method. 



F. two-species driven lattice gas 

In all the above examples, each site has only two local 
states (spin up or down) . In the two-species driven lattice 
gas model plj | , motivated in part by multi-ionic conduc- 
tors and traffic flow problems, there are three possibili- 
ties: as a hole, or either of two types of particles. The 



two types of particles are driven in opposite directions as 
if they were oppositely charged and driven by an electric 
field, with local particle densities denoted by p + and p- . 
Due to the extra local state, it is not easy to write down 
the correct set of equations by symmetry and intuition 
alone. One way to proceed |21| is to express the entropy 
in terms of p + and p_ , and obtain the diffusion terms by 
functional differentiations. The driving terms can then 
be added to the kinetic equations by generalizing that for 
the one-species model. 

Another way is to apply the current method |2^] . The 
equations derived are the same as in pT| , with the advan- 
tage of tractable microscopic origins in each coefficients. 
Hence, this model further testifies the usefulness of the 
present approach when symmetry and intuition are not 
very helpful. 

We end this section with a final remark. Since the 
method begins with factorization of joint probabilities, 
the ensuing equation is deterministic, all information 
about correlations seem to have lost. The situation can 
be remedied, however, by introducing a noise term to re- 
store a probabilistic description. Correlations can then 
be computed by averaging over the noise by means of 
standard field-theoretic techniques p3| . For equilibrium 
systems, the noise can be fixed by requirements such as 
the fluctuation-dissipation theorem. For noncquilibrium 
systems, there is no general rule. One usually has to 
extrapolate by analogy to equilibrium. 



IV. CONCLUSION 

We have presented in details a very simple and straight 
forward method to derive the deterministic kinetic equa- 
tions from known microscopic dynamics for stochastic, 
interacting systems. The method has a mean-field fla- 
vor. It preserves the underlying symmetries of the dy- 
namics and is in line with the spirit of coarse graining. 
The resulting equations are in good agreement with other 
either more or less rigorous approaches, as demonstrated 
explicitly via several examples. Hence, despite the ap- 
proximate and heuristic nature of our approach, it proves 
to be a useful and convenient means to obtain a correct 
continuum theory, especially (i) when other more rigor- 
ous approaches do not apply; (ii) when symmetries of 
the system is not intuitively obvious; and (iii) when mi- 
croscopic dependences of the continuum parameters are 
wanted. 

Acknowledgments: This work is supported by the Na- 
tional Science Council of R.O.C. under grant number 
NSC89-2112-M-001-015, and a Main-Theme Grant of the 
Academia Sinica. 
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APPENDIX A: NOISE CORRELATION 

Since the present method only gives the deterministic 
part of the kinetic equation, the noise term has to be con- 
sidered separately. Here we follow the common practice 
to assume that the noise £ is Gaussianly distributed and 
correlated over negligible ranges. Then the only question 
is to determine its correlation D in (Q). 

There are two ways to do it. The first makes use of 
the correspondence between the Langevin equation 



N.G. van Kampen, Stochastic Processes in 
Chemistry, (North-Holland, Amsterdam, 



dt 



-r 



sn 



C 



(C(f,t)t(f,1?)) = 2D6(r-7)S(t-t>), 
and the Fokker-Planck equation 



(Al) 
(A2) 



-m = -J dx rA W' (A3) 

which is a continuity equation. A stationary solution of 
the Fokker-Planck equation is obtained by setting zero 
the probability current, i.e., (•••) = 0, which gives V oc 
g-r-H/D g mce t ne free energy F[h] in the presence of an 
external field h is of the form T[h] = JF[0] — h(f>, it differs 
from Ti by a factor of /i. Hence, by matching g-^MA-BT 
and e~ rn ^ h ^ D , we deduce that 



D = k B TfiT. 



In passing, it is worth noting that the kinetic coefficient 
defined in 



at 



is A = D/ksT: the Einstein relation. 

An alternative way to determine D is to use the 
fluctuation-dissipation theorem, which in momentum- 
frequency space takes the form 



2k B T 



lmx(k,Lu) = G(k,u). 



Although neither the susceptibility \ nor the- two-point 
correlation function G can be calculated in closed form 
for general Ti, ( |A6| ) holds order by order so that we only 
need to consider the Gaussian model in the case of g — 
0. Thus, by Fourier transforms, we obtain x(fc,c<j) = 
Tfi/[-iuj+T(k 2 +r)} andG(Jfe,w) = 2_D/[o; 2 +r 2 (fc 2 +r) 2 ], 
which by virtue of ( |A6|) also gives ( A4). 
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FIG. 1. Quadratic coefficient in the transverse direction 
r x plotted vs. temperature. Its zero locates the critical tem- 
perature T C GL . 




FIG. 2. (a) Trends of r y vs. r x as T is varied across 
at fixed E. From bottom to top: E = 0, 2, 4, 6, 10, 14, 20, 
and 50. (b) Intercept r y (r x = 0) in (a) plotted vs. E. 
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FIG. 3. The coefficient of the leading nonlinearity, £ vs. 
the microscopic drive E at different temperatures. From top 
to bottom: k B T / J = 1, 3, 3.86143(= T C GL ), 5, 10, 20 and 50. 
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